package MeshlessShapeFunctions;

import HashMatrix2.HashMatrix2;
import HashMatrix2.LUFactorization;
import Integration.IFunction2D;
import java.util.ArrayList;

/*
 * To change this template, choose Tools | Templates
 * and open the template in the editor.
 */

/**
 *
 * @author mark
 */
public class JMeshlessInterpolationFunction implements IFunction2D{
    protected int PolinomialDegree;
    protected double xii;
    protected double yii;
    protected double rho;
    protected IKernelFunction func;


    ArrayList<Double> BaseVector;
    LUFactorization lu;
    double[] alpha;
    boolean key;

    protected HashMatrix2 MBase;
    protected HashMatrix2 A;
    protected HashMatrix2 Ax;
    protected HashMatrix2 Ay;
    protected double[] b;
    protected ArrayList<Particle> NodeList;


    public JMeshlessInterpolationFunction(int PolinomialDegree, double xi, double yi, double rho, IKernelFunction func) {
        this.PolinomialDegree = PolinomialDegree;
        this.xii = xi;
        this.yii = yi;
        this.rho = rho;
        this.func = func;
    }

    public void setNodeList(ArrayList<Particle> NodeList) {
        this.NodeList = NodeList;

        for (int k = 0; k < NodeList.size(); k++) {
            Particle particle = NodeList.get(k);
            ArrayList<Double> BaseVectorTmp = getBaseVector(particle.x , particle.y,particle.x , particle.y);
            int n = BaseVectorTmp.size();
            MBase = new HashMatrix2(HashMatrix2.ARRAY_TYPE, n);

            for (int i = 0; i < n; i++) {
                double val = BaseVectorTmp.get(i);
                for (int j = 0; j < n; j++) {
                    MBase.addValue(i, j, val*BaseVectorTmp.get(j));
                }
            }
        }
    }

    public ArrayList<Double> getBaseVector(double x,double y,double xi, double yi){
        ArrayList<Double> L = new ArrayList<Double>();
        for (int i = 0; i <= this.PolinomialDegree; i++) {
            for (int j = 0; i+j <= this.PolinomialDegree; j++) {
                L.add(myPow(xi, i)*myPow(yi, j));
            }
        }
        return L;
    }

    private double myPow(double x,double i){
        if(i>=0.0){
            return Math.pow(x, i);
        }else{
            return 0.0;
        }
    }

    public ArrayList<Double> getBaseVectorDerivatyX(double x,double y,double xi, double yi){
        ArrayList<Double> L = new ArrayList<Double>();
        for (int i = 0; i <= this.PolinomialDegree; i++) {
            for (int j = 0; i+j <= this.PolinomialDegree; j++) {
                L.add((i)*myPow(xi, i-1)*myPow(yi, j));
            }
        }
        return L;
    }

    public ArrayList<Double> getBaseVectorDerivatyY(double x,double y,double xi, double yi){
        ArrayList<Double> L = new ArrayList<Double>();
        for (int i = 0; i <= this.PolinomialDegree; i++) {
            for (int j = 0; i+j <= this.PolinomialDegree; j++) {
                L.add((j)*myPow(xi, i)*myPow(yi, j-1));
            }
        }

        return L;
    }

    public double calculate(double x, double y){
        BaseVector = getBaseVector(x, y, xii, yii);
        calculateAlpha(x, y);
        double returnValue = 0.0;

        if(key){
            lu = new LUFactorization(A);
            alpha = lu.calculate(b);

            BaseVector = getBaseVector(x, y, x, y);

            for (int i = 0; i < b.length; i++) {
                returnValue+=BaseVector.get(i)*alpha[i];
            }
        }
        alpha = null;
        BaseVector = null;
        lu=null;
        return returnValue;
    }

    private void calculateAlpha(double x,double y){
        BaseVector = getBaseVector(x, y, xii, yii);
        int n = BaseVector.size();
        b = new double[n];
        A = new HashMatrix2(HashMatrix2.ARRAY_TYPE, n);
        for (int i = 0; i < NodeList.size(); i++) {
            Particle particle = NodeList.get(i);
            calculateAlphaInside(x,y,particle.x, particle.y);
        }
        key = false;
        for (int i = 0; i < b.length; i++) {
            b[i] = BaseVector.get(i)*func.evaluateSque((xii-x)/rho, (yii-y)/rho);//////////////
            if(b[i]!=0.0){key=true;}
        }
    }

    private void calculateAlphaInside(double x,double y,double xi, double yi){
        ArrayList<Double> BaseVectorTmp = getBaseVector(xi, yi, xi, yi);
        int n = BaseVectorTmp.size();
        for (int i = 0; i < n; i++) {
            double val = BaseVectorTmp.get(i)*func.evaluateSque((xi-x)/rho, (yi-y)/rho);/////////////////////
            for (int j = 0; j < n; j++) {
                A.addValue(i, j, val*BaseVectorTmp.get(j));
            }
        }
    }

    private void calculateAlfaX(double val){
        int n = MBase.getN();
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < n; j++) {
                Ax.setValue(i, j, MBase.getValue(i, j)*val);
            }
        }
    }

    public double calculateDerivatyX(double x, double y) {
        ArrayList<Double> BaseVectorTmp = getBaseVector(x, y, xii, yii);
        int n = BaseVectorTmp.size();
        b = new double[n];
        Ax = new HashMatrix2(HashMatrix2.ARRAY_TYPE, n);
        double val = 0.0;
        for (int i = 0; i < NodeList.size(); i++) {
            Particle particle = NodeList.get(i);
            val+=func.evaluatePrimeNormX((particle.x-x)/rho, (particle.y-y)/rho);/////////////////////////
        }

        calculateAlfaX(val);

        //if(lu==null){
            //if(A==null){
                this.calculateAlpha(x, y);
            //}
            lu = new LUFactorization(A);
            alpha = lu.calculate(b);
        //}

        key = false;

        ArrayList<Double> BaseVectorTmp1 = this.getBaseVectorDerivatyX(x, y, xii, yii);

        for (int i = 0; i < b.length; i++) {
            double tmp = 0.0;
            for (int j = 0; j < n; j++) {
                tmp+=alpha[j]*Ax.getValue(i, j);
            }
            b[i] = BaseVectorTmp1.get(i) + tmp;
            if(b[i]!=0.0){key=true;}
        }


        double returnValue = 0.0;

        if(key){
            //lu = new LUFactorization(Ax);
            double[] res = lu.calculate(b);

            for (int i = 0; i < b.length; i++) {
                returnValue+=BaseVectorTmp.get(i)*(res[i]*func.evaluateNorm((xii-x)/rho, (yii-y)/rho) + alpha[i]*func.evaluatePrimeNormX((xii-x)/rho, (yii-y)/rho));/////////////////////////
            }
        }
        if(String.valueOf(returnValue).equals("Nan")){return 0.0;}
        return returnValue;
    }

    public double calculateDerivatyY(double x, double y) {
        throw new UnsupportedOperationException("Not supported yet.");
    }

    public double calculateDiffuseDerivatyX(double x,double y){
        BaseVector = getBaseVector(x, y, xii, yii);
        calculateAlpha(x, y);
        double returnValue = 0.0;

        if(key){
            lu = new LUFactorization(A);
            alpha = lu.calculate(b);

            BaseVector = getBaseVectorDerivatyX(x, y, x, y);

            returnValue = 0.0;
            for (int i = 0; i < alpha.length; i++) {
                returnValue+=BaseVector.get(i)*alpha[i];
            }
        }

        alpha = null;
        BaseVector = null;
        lu=null;

        return returnValue;
    }

    public double calculateDiffuseDerivatyY(double x,double y){
        BaseVector = getBaseVector(x, y, xii, yii);
        calculateAlpha(x, y);
        double returnValue = 0.0;

        if(key){
            lu = new LUFactorization(A);
            alpha = lu.calculate(b);

            BaseVector = getBaseVectorDerivatyY(x, xii, y, yii);

            returnValue = 0.0;
            for (int i = 0; i < alpha.length; i++) {
                //System.out.println(BaseVector.get(i)+"\t\t\t"+alpha[i]);
                returnValue+=BaseVector.get(i)*alpha[i];
            }
            //System.out.println("================================================================");
        }

        alpha = null;
        BaseVector = null;
        lu=null;

        return returnValue;
    }

    public double getRho() {
        return rho;
    }

    public double getXii() {
        return xii;
    }

    public double getYii() {
        return yii;
    }

    

    public static void main(String[] args){
        ArrayList<Particle> nodes = new ArrayList<Particle>();
//        ArrayList<JMeshlessInterpolationFunction> listBaseFunction = new ArrayList<JMeshlessInterpolationFunction>();
        //параметри області:
        double ax = 0, bx = 1,
               ay = 0, by = 1;

        int n = 5, m = 5;

        double  hx = (bx-ax)/(double)n,
                hy = (by-ay)/(double)m;

        //вузли+функції:
        for (int i = 0; i <= n; i++) {
            for (int j = 0; j <= m; j++) {
                Particle p = new Particle(ax+i*hx, ay+j*hy, i*n+j);
                nodes.add(p);
            }
        }
        JMeshlessInterpolationFunction f = new JMeshlessInterpolationFunction(3, 0.5, 0.5, 0.65, new JGaussian(0.5));
        f.setNodeList(nodes);

        nodes = new ArrayList<Particle>();
        n = 25; m = 25;

                hx = (bx-ax)/(double)n;
                hy = (by-ay)/(double)m;

        //вузли+функції:
        for (int i = 0; i <= n; i++) {
            for (int j = 0; j <= m; j++) {
                Particle p = new Particle(ax+i*hx, ay+j*hy, i*n+j);
                nodes.add(p);
            }
        }
        for (int i = 0; i < nodes.size(); i++) {
            Particle particle = nodes.get(i);
            System.out.println(particle.x+"\t"+particle.y+"\t"+f.calculateDiffuseDerivatyY(particle.x, particle.y));
        }

    }

}

//package MeshlessShapeFunctions;
//
//import HashMatrix2.HashMatrix2;
//import HashMatrix2.LUFactorization;
//import java.util.ArrayList;
//
///*
// * To change this template, choose Tools | Templates
// * and open the template in the editor.
// */
//
///**
// *
// * @author mark
// */
//public class JMeshlessInterpolationFunction {
//    protected int PolinomialDegree;
//    protected double xii;
//    protected double yii;
//    protected double rho;
//    protected IKernelFunction func;
//    protected HashMatrix2 A;
//    protected double[] b;
//    protected ArrayList<Particle> NodeList;
//
//
//    public JMeshlessInterpolationFunction(int PolinomialDegree, double xi, double yi, double rho, IKernelFunction func) {
//        this.PolinomialDegree = PolinomialDegree;
//        this.xii = xi;
//        this.yii = yi;
//        this.rho = rho;
//        this.func = func;
//    }
//
//    public void setNodeList(ArrayList<Particle> NodeList) {
//        this.NodeList = NodeList;
//    }
//
//    public double getRho() {
//        return rho;
//    }
//
//    public void setRho(double rho) {
//        this.rho = rho;
//    }
//
//    public double getXii() {
//        return xii;
//    }
//
//    public void setXii(double xii) {
//        this.xii = xii;
//    }
//
//    public double getYii() {
//        return yii;
//    }
//
//    public void setYii(double yii) {
//        this.yii = yii;
//    }
//
//    public ArrayList<Double> getBaseVector(double x,double y,double xi, double yi){
//        ArrayList<Double> L = new ArrayList<Double>();
//        for (int i = 0; i <= this.PolinomialDegree; i++) {
//            for (int j = 0; i+j <= this.PolinomialDegree; j++) {
//                L.add(Math.pow((xi-x)/rho, i)*Math.pow((yi-y)/rho, j));
//            }
//        }
//
//        return L;
//    }
//
//    public ArrayList<Double> getBaseVectorDerivatyX(double x,double y,double xi, double yi){
//        ArrayList<Double> L = new ArrayList<Double>();
//        for (int i = 0; i <= this.PolinomialDegree; i++) {
//            for (int j = 0; i+j <= this.PolinomialDegree; j++) {
//                L.add(Math.pow((xi-x)/rho, i)*Math.pow((yi-y)/rho, j));
//            }
//        }
//
//        return L;
//    }
//
//        public ArrayList<Double> getBaseVectorDerivatyY(double x,double y,double xi, double yi){
//        ArrayList<Double> L = new ArrayList<Double>();
//        for (int i = 0; i <= this.PolinomialDegree; i++) {
//            for (int j = 0; i+j <= this.PolinomialDegree; j++) {
//                L.add(Math.pow((xi-x)/rho, i)*Math.pow((yi-y)/rho, j));
//            }
//        }
//
//        return L;
//    }
//
//    public double calculate(double x, double y){
//        if(((x-xii)*(x-xii)+(y-yii)*(y-yii))>rho*rho){return 0;}
//        ArrayList<Double> BaseVector = getBaseVector(x, y,xii,yii);
//        int n = BaseVector.size();
//        b = new double[n];
//        A = new HashMatrix2(HashMatrix2.ARRAY_TYPE, n);
//        for (int i = 0; i < NodeList.size(); i++) {
//            Particle particle = NodeList.get(i);
//            calculateAlfa(x,y,particle.x,particle.y);
//        }
//
//        boolean key=false;
//        for (int i = 0; i < b.length; i++) {
//            b[i] = BaseVector.get(i);
//            if(b[i]!=0.0){key = true;}
//        }
//
//        double returnValue = 0.0;
//        if(key){
//            LUFactorization lu = new LUFactorization(A);
//            double[] res = lu.calculate(b);
//
//            BaseVector = getBaseVector(x, x,y,y);
//            for (int i = 0; i < b.length; i++) {
//                returnValue+=BaseVector.get(i)*res[i]*func.evaluateSque((xii-x)/rho,(yii-y)/rho);
//            }
//        }
//        return returnValue;
//    }
//
//
//    private void calculateAlfa(double x,double y,double xi, double yi){
//        ArrayList<Double> BaseVector = getBaseVector(x, y, xi, yi);
//        int n = BaseVector.size();
//        //b = new double[n];
//        for (int i = 0; i < n; i++) {
//            double val = BaseVector.get(i)*func.evaluateSque((xi-x)/rho, (yi-y)/rho);
//            for (int j = 0; j < n; j++) {
//                A.addValue(i, j, val*BaseVector.get(j));
//            }
//        }
//    }
//    public double calculateDiffuseDerivatyX(double x, double y){
//        ArrayList<Double> BaseVector = getBaseVector(x,y, xii,yii);
//        int n = BaseVector.size();
//        b = new double[n];
//        A = new HashMatrix2(HashMatrix2.ARRAY_TYPE, n);
//        for (int i = 0; i < NodeList.size(); i++) {
//            Particle particle = NodeList.get(i);
//            calculateAlfa(x,y,particle.x,particle.y);
//        }
//
//        boolean key=false;
//        for (int i = 0; i < b.length; i++) {
//            b[i] = BaseVector.get(i)*func.evaluateSque((xii-x)/rho,(yii-y)/rho);
//            if(b[i]!=0.0){key = true;}
//        }
//
//        double returnValue = 0.0;
//        if(key){
//            LUFactorization lu = new LUFactorization(A);
//            double[] res = lu.calculate(b);
//
//            BaseVector = getBaseVectorDerivatyX(x, y,x,y);
//
//            for (int i = 0; i < b.length; i++) {
//                returnValue+=BaseVector.get(i)*res[i];
//            }
//        }
//        return returnValue;
//    }
//
//    public double calculateDiffuseDerivatyY(double x, double y){
//        ArrayList<Double> BaseVector = getBaseVector(x,y, xii,yii);
//        int n = BaseVector.size();
//        b = new double[n];
//        A = new HashMatrix2(HashMatrix2.ARRAY_TYPE, n);
//        for (int i = 0; i < NodeList.size(); i++) {
//            Particle particle = NodeList.get(i);
//            calculateAlfa(x,y,particle.x,particle.y);
//        }
//
//        boolean key=false;
//        for (int i = 0; i < b.length; i++) {
//            b[i] = BaseVector.get(i)*func.evaluateSque((xii-x)/rho,(yii-y)/rho);
//            if(b[i]!=0.0){key = true;}
//        }
//
//        double returnValue = 0.0;
//        if(key){
//            LUFactorization lu = new LUFactorization(A);
//            double[] res = lu.calculate(b);
//
//            BaseVector = getBaseVectorDerivatyX(x, y,x,y);
//
//            for (int i = 0; i < b.length; i++) {
//                returnValue+=BaseVector.get(i)*res[i];
//            }
//        }
//        return returnValue;
//    }
//
//    public double[] calculateAll(double x, double y){
//        ArrayList<Double> BaseVector = getBaseVector(x,y, xii,yii);
//        int n = BaseVector.size();
//        b = new double[n];
//        A = new HashMatrix2(HashMatrix2.ARRAY_TYPE, n);
//        for (int i = 0; i < NodeList.size(); i++) {
//            Particle particle = NodeList.get(i);
//            calculateAlfa(x,y,particle.x,particle.y);
//        }
//
//        boolean key=false;
//        for (int i = 0; i < b.length; i++) {
//            b[i] = BaseVector.get(i)*func.evaluateSque((xii-x)/rho,(yii-y)/rho);
//            if(b[i]!=0.0){key = true;}
//        }
//
//        double[] returnV = new double[3];
//        double returnValue = 0.0;
//        if(key){
//            LUFactorization lu = new LUFactorization(A);
//            double[] res = lu.calculate(b);
//
//            BaseVector = getBaseVector(x, y,x,y);
//
//            for (int i = 0; i < b.length; i++) {
//                returnValue+=BaseVector.get(i)*res[i];
//            }
//            returnV[0] = returnValue;
//
//            BaseVector = getBaseVectorDerivatyX(x, y,x,y);
//
//            for (int i = 0; i < b.length; i++) {
//                returnValue+=BaseVector.get(i)*res[i];
//            }
//            returnV[1] = returnValue;
//
//            BaseVector = getBaseVectorDerivatyY(x, y,x,y);
//
//            for (int i = 0; i < b.length; i++) {
//                returnValue+=BaseVector.get(i)*res[i];
//            }
//            returnV[2] = returnValue;
//        }
//
//        return returnV;
//    }
//
//    public static void main(String[] args){
//        ArrayList<Particle> nodes = new ArrayList<Particle>();
//        ArrayList<JMeshlessInterpolationFunction> listBaseFunction = new ArrayList<JMeshlessInterpolationFunction>();
//        //параметри області:
//        double ax = 0, bx = 1,
//               ay = 0, by = 1;
//
//        int n = 25, m = 25;
//
//        double  hx = (bx-ax)/(double)n,
//                hy = (by-ay)/(double)m;
//
//        //вузли+функції:
//        for (int i = 0; i < n; i++) {
//            for (int j = 0; j < m; j++) {
//                Particle p = new Particle(ax+i*hx, ay+j*hy, i*n+j);
//                nodes.add(p);
//            }
//        }
//        JMeshlessInterpolationFunction f = new JMeshlessInterpolationFunction(2, 0.5, 0.5, 0.5, new JGaussian());
//        f.setNodeList(nodes);
//        for (int i = 0; i < nodes.size(); i++) {
//            Particle particle = nodes.get(i);
//            System.out.println(particle.x+"\t"+particle.y+"\t"+f.calculate(particle.x, particle.y));
//        }
//
//    }
//
//}
